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Nuclear matter and compact neutron stars are studied in the framework of the non-linear derivative (NLD) 
model which accounts for the momentum dependence of relativistic mean-fields. The generalized form of 
the energy-momentum tensor is derived which allows to consider different forms of the regulator functions 
in the NLD Lagrangian. The thermodynamic consistency of the NLD model is demonstrated for arbitrary 
choice of the regulator functions. The NLD approach describes the bulk properties of the nuclear matter and 
compares well with microscopic calculations and Dirac phenomenology. We further study the high density 
domain of the nuclear equation of state (EoS) relevant for the matter in /3-equilibrium inside neutron stars. It 
is shown that the low density constraints imposed on the nuclear EoS and by the momentum dependence of the 
Schrodinger-equivalent optical potential lead to a maximum mass of the neutron stars around M ~ 2Mq which 
accommodates the observed mass of the J1614-2230 millisecond radio pulsar. 
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I. INTRODUCTION 

Relativistic hadrodynamics (RHD) of interacting nucleons 
and mesons provide a simple and successful tool for the theo- 
retical description of different nuclear systems such as nuclear 
matter, finite nuclei, heavy-ion collisions and compact neutron 
stars [1]. Starting from the pioneering work of Duerr [2], sim- 
ple RHD Lagrangians have been introduced [3, 4] and since 
then many different extensions of RHD approach, which rely 
on relativistic mean-field (RMF) approximation, have been 
developed. They describe the saturation mechanism in nu- 
clear matter and generate a natural mechanism for the strong 
spin-orbit force in nuclei. An energy dependence of the 
Schrodinger-equivalent optical potential [5, 6] is thereby in- 
cluded as a consequence of a relativistic description. How- 
ever, when using the standard RHD Lagrangian in RMF ap- 
proximation, the nucleon selfenergies become simple func- 
tions of density only, and do not depend explicitly on momen- 
tum of the nucleon. As a consequence a linear energy de- 
pendence of the Schrodinger-equivalent optical potential with 
a divergent behavior at high energies arises [7]. This well- 
known feature contradicts Dirac phenomenology [5, 6, 8]. 

To solve this issue one may go beyond the mean-field ap- 
proximation in a quantum field theoretical framework by a 
systematic diagrammatic expansion of nucleon selfenergies. 
For instance, in Dirac-Brueckner-Hartree-Fock (DBHF) [9- 
11] calculations the nucleon selfenergies indeed depend on 
both the density and single particle momentum. They repro- 
duce the empirical saturation point of nuclear matter as well as 
the energy dependence of the optical potential at low energies. 
However, the DBHF approach has its apparent limitations at 
high energies and densities relevant, for instance, in heavy-ion 
collisions where its application within transport theory turns 
out to be intricate [12, 13]. Also the thermodynamic consis- 
tency of the DBHF calculations is not obvious [14]. 
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As an alternative approach to ab-initio DBHF calculations 
for the nuclear many-body systems a phenomenological treat- 
ment of the problem in the spirit of the RMF approximation 
is still considered as a powerful tool. However, the simple 
Lagrangian of RHD [3, 4] has to be further modified for a 
quantitative description of static nuclear systems such as nu- 
clear matter and/or finite nuclei. Therefore, it is mandatory 
to introduce new terms, e.g., including non-linear self inter- 
actions of the scalar [15] and/or vector [16] meson fields, or 
to modify existing contributions in the Lagrangian, e.g., intro- 
ducing density dependent meson-nucleon couplings [17-19]. 
The model parameters have to be then fitted to properties of 
nuclear matter and/or atomic nuclei, since, they cannot be de- 
rived in a simple manner from a microscopic description. 

The momentum dependence of in-medium interactions be- 
comes particularly important in description of nuclear colli- 
sion dynamics such as heavy-ion collisions. Indeed, analyses 
of proton-nucleus scattering data [5, 6] show that the proton- 
nucleus optical potential starts to level off already at incident 
energies of about 300 Me V. Thus, other RMF approaches have 
been developed by including additional non-local contribu- 
tions, i.e., by introducing Fock-terms, on the level of the RMF 
selfenergies leading to a density and momentum dependent 
interactions [7]. However, such a treatment is not covariant 
and also its numerical realization in actual transport calcula- 
tions is rather difficult [7]. Another approach has been pro- 
posed in [20] and more recently in [21, 22] by introducing 
higher order derivative couplings in the Lagrangian of RHD. 
In Ref. [20] such gradient terms have been studied with the 
conclusion of a softening of the nuclear EoS. In another study 
of Ref. [22] both the density dependence of the nuclear EoS 
and the energy dependence of the optical potential have been 
investigated. While the modified interactions of meson fields 
with nucleons explain the empirical energy dependence of the 
optical potential, a stiff EoS at high densities results from an 
introduction of an explicit density dependence of the nucleon- 
meson couplings with additional parameters. The impact of 
momentum dependent RMF models on nuclear matter bulk 
properties and particularly on the high density domain of EoS 
relevant for neutron stars is presently less understood. 
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The purpose of the present work is to develop a relativistic 
and thermodynamically consistent RMF model, which pro- 
vides the correct momentum dependence of the nucleon self- 
energies and agrees well with available empirical information 
on nuclear matter ground state, in a self consistent Lagrangian 
framework. Some steps in this direction have been already 
done in Refs. [23-25] where the concept of non-linear deriva- 
tive meson-nucleon Lagrangian has been introduced. How- 
ever, the calculations of Refs. [23-25] were based on a partic- 
ular exponential form of the regulators in the RHD Lagrangian 
and a detailed study of nuclear matter ground state proper- 
ties has not been done. In the present work the generalized 
form of the energy-momentum tensor in the NLD model is 
derived and allows to consider different regulator functions in 
the Lagrangian. The thermodynamic consistency of the NLD 
model is demonstrated for arbitrary choice of the regulators. 
A thorough study of the properties of nuclear matter around 
saturation density is further performed. The model describes 
the bulk properties of the nuclear matter and compares well 
with microscopic calculations and Dirac phenomenology. We 
also investigate the high density region of the NLD EoS rel- 
evant for the neutron stars. It is found that the low density 
constraints imposed on the nuclear matter EoS and by the 
momentum dependence of the Schrodinger-equivalent opti- 
cal potential lead to a maximum mass of the neutron stars 
around M ~ 2M Q . It is demonstrated that the high den- 
sity pressure-density diagram as extracted from astrophysical 
measurements [26, 27] can be well described with nucleonic 
degrees of freedom only. 



II. FIELD THEORY WITH HIGHER DERIVATIVES 

The non-linear derivative (NLD) model is based on a field- 
theoretical formalism which accounts for the higher-order 
derivative interactions in the RHD Lagrangian. As a conse- 
quence, the conventional RHD mean-field theory based on 
minimal interaction Lagrangians has to be extended to the 
case of higher-order non-linear derivative functionals. 

For that purpose we consider the most general structure of 
a Lagrangian density C with higher-order field derivatives, i.e. 

C((f r (x), d ai tp r (x), d aia2 (p r (x),---,d ai ... an (p r (x)) , (1) 

where it is supposed that C has continuous derivatives up to 
order n with respect to all its arguments, that is 
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dx ai dx 



— tpr(x) = 8 ai ■ ■ ■ d an (fir(x) 



where on is a four index and x denotes the coordinates in 
Minkowski space. The order n can be a finite number or 
n — > oo. The subscript r denotes different fields, for instance, 
in the case of the spinor fields one would have <pi = $ and 

The derivation of the generalized Euler-Lagrange equations 
of motion follows from the variation principle for the action 
S = J d 4 xC(x) with the Lagrangian of Eq. (1), where one 
considers tp r , d ai tp r , d ai d a2 <p r , d ai - ■ ■d oln ip r as inde- 
pendent generalized coordinates. The Euler-Lagrange equa- 



tions are obtained from principle of least action 

5S = , (2) 

where SS is given by 

SS = J d 4 x5C (ip r , d ai (p r , d aia2 ip r , ■ ■ ■, d ai ... an ip r ) (3) 
and is obtained by the variation of the generalized coordinates 



<Pr 
d ai tp r 



lfr + 5 ip r , 
d ai lfi r +5d ai lfi r , 



d ai . 



% (p r + Sd ai . 



(4) 



with vanishing contributions on the surface of the integration 
volume as the boundary condition. The variation of the La- 
grangian density with respect to all degrees of freedom reads 
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dC . dC 

dip r o(d ai ip r 
dC 



■d ai S(fi r + 



dC 



d{d aia2 ip r ) 



d{d ai ... an ip r ) 



d aia2 5ip r 



(5) 



As a next step one inserts Eq. (5) into Eq. (3) and then per- 
forms successively partial integrations, e.g., one partial inte- 
gration for the second term in Eq. (5), two partial integrations 
for the third term in Eq. (5), and n partial integrations for the 
last term. This procedure results in to the following integrand 
in Eq. (3) 
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d£ L _ d dC 

dlfr ^ d(d ai ip r ) 

+ ■ •• + (-)" d ai ... c 



dC 



dC 



d{d aia2 (f r ) 



d(d ai ... an Ifr) 



S(f r (6) 



up to 4-divergence terms, which by Gauss law do not con- 
tribute to the action in Eq. (3). Thus, one arrives to the fol- 
lowing generalized Euler-Lagrange equation 



dC 

dcp r 



dC 



i=l 



d(d ai ... ai <p r ) 



(7) 



The Noether theorem follows from invariance principles of 
the Lagrangian density, Eq. (1), with respect to infinitesimal 
variations of the generalized coordinates and their argument 
x^ (see for notations Appendix A). As further shown in Ap- 
pendix B, the requirement of invariance of the Lagrangian 
density, Eq. (1), with respect to global phase transformations 



(fr{x) > (f' r (x) = e te (p r {x) 



(8) 



leads to a continuity equation <9 M J M = for a conserved 
Noether current J M . The latter is given by the following ex- 
pression 
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= -l [K$<p r + K^d^Pr + K^d^ifr + ■■■ + d ai ...„„ <p r ] . (9) 

In fact, for n — > oo the Noether current consists of an infinite sequence of tensors with increasing rank order. Furthermore, each 
of the different tensors /C£," Tl<T2 "' in Eq. (9) contains again infinite series terms of higher-order derivatives with respect to the 
Lagrangian density. They are given by the following expressions 

n i-1 fjr 

n i—1 r^p 

^r=E(-) i+i n^^7) = 

n i — 1 ^ „ 



The derivation of the energy-momentum tensor proceeds in a similar way, see Appendix B. Now the field arguments are 
transformed, but not the fields them self. In particular, invariance of the Lagrangian density (1) with respect to a constant 
displacement 5^ of the coordinates x^ 

— > x'^ = Xf, + 8^ , (11) 
implies a continuity equation d^T^ = for the energy-momentum tensor T^ v which takes the following form 

T» v = K»d v <p r + KUTKiVt + 1C^d»^ip r + ■■■ + K^-^d" ai ... a ^ r - gTC . (12) 

The 00-component of the energy-momentum tensor describes the energy density and the spatial diagonal components are related 
to the pressure density. These equations form a background for the construction and application of the NLD formalism presented 
in the proceeding sections. They will further provide a thermodynamically consistent framework for the calculation of the EoS 
in mean field approximation in terms of energy and pressure densities. 



III. THE NON-LINEAR DERIVATIVE MODEL 

In this section we introduce the non-linear derivative (NLD) model and derive the equations of motion for the relevant degrees 
of freedom. The NLD approach is essentially based on the Lagrangian density of RHD [2^], which is given by 
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dj d"S + & 



(13) 



where * = (* p , *„) T denotes the nucleon spinor field in the 
Lagrangian density of a Dirac-type. In a spirit of RHD, the 
interactions between the nucleon fields are described by the 
exchange of meson fields. These are the scalar a and vector 
mesons in the isoscalar channel, as well as the scalar 8 
and vector p M mesons in the isovector channel. Their corre- 
sponding Lagrangian densities are of the Klein-Gordon and 
Proca types, respectively. The term U{a) = \b<7^ + jca 4 



contains the usual selfinteractions of the a meson. The no- 
tations for the masses of fields in Eq. (13) are obvious. The 
field strength tensors are defined as = d^u) v — d v uj^, 
G ^ = d^p v — d v ' p^ for the isoscalar and isovector fields, 
respectively. 



In conventional RHD approaches the interaction La- 
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grangian C lnt is given by [3, 4] 

r per I ru> \ rP i pS 

>-int '-int i '-'int ~ '-int < l -int ) 

where 



(14) 



P _ _9p 



(20) 



(21) 



cLt = -g P ^r^p t 



(15) 
(16) 
(17) 
(18) 



and C 



cr,oj,p,8 



contains the meson-nucleon interactions with 



coupling strengths g a .u. P .s an d T denotes the isospin Pauli op- 
erator. 

In the NLD model the momentum dependence of fields is 
realized by the introduction of non-linear derivative operators 
in the interaction Lagrangian of conventional RHD. These ad- 
ditional operators regulate the high momentum components 
of the RMF fields in the interaction vertices and can be in- 
terpreted as cut-off form factors. This is in spirit of boson- 
exchange models where the phenomenological cut-off is an 
indispensable part of any microscopic description of meson- 
nucleon interaction [28, 29]. In the RMF (Hartree) approx- 
imation to RHD only bare Lorentz structures corresponding 
to the point-like meson-nucleon interactions are taken into 
account and the high momentum components of fields are 
not suppressed due to the missing nucleon finite size effect. 
The NLD model attempts to account for the suppression of 
the high momentum part of the nucleon field in the meson- 
nucleon interaction on a field-theoretical level. 

The NLD interaction Lagrangians contain the conventional 
meson-nucleon RHD structures, however, they are extended 
by the inclusion of non-linear derivative operators into the 
meson-nucleon vertices. The NLD interaction Lagrangians 
followed here read 



pa ao_ 

int r) 



(19) 



f s = 9s_ 



# W*<5 + <f*f 



(22) 



As one can see, the only difference with respect to the con- 
ventional RHD interaction Lagrangian is the presence of ad- 
ditional operators T> , V which serve to regulate the high mo- 
mentum component of the nucleon field. The hermiticity of 
the Lagrangian demands The operator functions 

(regulators) are assumed to be generic functions of 

partial derivative operator and supposed to act on the nucleon 
spinors "J and "J, respectively. Furthermore, these regulators 
are assumed to be smooth functions. Therefore, the formal 
Taylor expansion of the operator functions in terms of partial 
derivatives generates an infinite series of higher-order deriva- 
tive terms 



V :=£>( £ 



=^1^0 \ ' (23) 



E 

j 



3=0 



(24) 



The expansion coefficients are given by the partial derivatives 
of V with respect to the operator arguments £ and^ around 
the origin. The operators are defined as £ ——( a id a , £ = 



* d a( a where the four vector £ p 



'/A contains the cut- 



off A and is an auxiliary vector. The functional form of 
the regulators is constructed such that in the limit A — > oo 
the following limit holds V( V) — > 1. Therefore, in the limit 
A — > oo the original RHD Lagrangians are recovered. 

In the most general case the NLD formalism can be ex- 
tended to the case of multiple variable regulators. In partic- 
ular, we can assume the non-linear operator to be a multi- 
variable non-linear function of higher-order partial deriva- 
tives, which are given by the following Taylor expansion 



^1,^2,^3,^4) 



E E 



i 1= i 2 =0 i 3 =0 i 4 =0 <5 £ £ 2® £ 3 3 ^ £ 4* 



E E E E 

ii=0 i2—0 i3—0 i4—0 



y ii y 12 \Ti3 
Si S 2 S3 S . 



-D\ _> _> _> _> 



S 1 S 2 S3 S 4 



(25) 



^\ f t t t t wo" < 26 > 



Then Eqs. (25) and (26) can be rearranged into the terms with way as before 
increasing order with respect to the partial derivatives, see for 

details Appendix C. The operators & are defined in a similar -> t «— >n 

ii = -Qid a , £i=i<9 Q Cf, (27) 
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with = v If /A (i = 1, 2, 3, 4) in this case. As we will show 
latter on, this representation allows to generate any desired 
form of the regulator function, i.e., momentum and/or energy 
dependent monopole, dipole etc. functions. 

The derivation of the equation of motion for the Dirac field 
follows the generalized Euler-Lagrange equations, Eq. (7), to 
the NLD-Lagrangian density using the Taylor form of the reg- 
ulators. This obviously will generate an infinite number of 
partial derivative terms in the equations of motions. However, 
as shown in detail in Appendix D these infinite series can be 
resummed (up to terms containing the derivatives of the me- 
son fields) to the following Dirac equation 



[^(id* - E") - (m - E,)] * = 



where the selfenergies E M and E s are given by 

E M = g^V + g p r ■ p^V +■■■ 
£ s = 9<y<^> + gsT - SV-\ . 



(28) 

(29) 
(30) 



Here both Lorentz-components of the selfenergy, S M and S s , 
show an explicit linear behavior with respect to the meson 
fields a, p M and 8 as in the standard RMF. However, they 
contain an additional dependence on regulator functions. 

The additional terms in Eqs. (29) and (30) containing the 
meson field derivatives are denoted by multiple dots. All these 
contributions can be also resummed. However, in the mean- 
field approximation to infinite nuclear matter, which will be 
discussed in the next section, these terms vanish. On the other 
hand, they will be needed in the description of finite systems, 



such as finite nuclei and heavy-ion collisions. Therefore, for 
simplicity we do not consider these terms here, and postpone 
the effect of these terms for future studies. 

The derivation of the meson field equations of motion is 
straightforward, since here one has to use the standard Euler- 
Lagrange equations 



— -d dc 

dif r " d(d a lfi r ) 



0, 



(31) 



where now r = a,u), p and 5. The following Proca and Klein- 
Gordon equations are obtained 



2 dU 1 
d a d a + m a a + — = -g, 

1 



dyGV +m 2 p p» 



d a d a 5 + m 2 J 



->9 P 



;96 



WW+WX^* , (32) 

* VY^ , (33) 

*^7"r 7^*1 , (34) 

f Pff + ffC$l . (35) 



Finally, we provide the general expressions for the Noether 
theorems within the NLD formalism. The evaluation of the 
conserved baryon current results from the application of the 
generalized expression for J M , Eq. (9), to the Lagrangian den- 
sity of the NLD model. As shown in detail in Appendix E, a 
systematic evaluation of the higher-order field derivatives of 
the NLD Lagrangian and the resummation procedure result in 



1 

+ 2-9P 



* t^7 Q * - *7 a 7^* 



■ 96 



(36) 



The new non-linear derivative operators in Eq. (36), fl^ := 
d i D/d{i ! d p ) and 7^ := d!$ /d(i~3 J, denote the deriva- 
tives of u and with respect to their operator argument i d M 

and id M (see Appendix E). The first term in Eq. (36) cor- 
responds to the standard expression of the RHD models and 



the additional contributions arise due to the additional higher- 
order field derivatives in the Noether theorem, Eq. (9). 

The energy-momentum tensor, T' 1!/ , is determined accord- 
ing to Eq. (12). The evaluation procedure, which is similar to 
that one for the Noether current, results in the following NLD 
expression for 



2 2 ' 



9a 



1 

"2 
1 

- 2 9 P 



v-^9u 



*7 Q 7^ i§ v * + # id v t?^7 a * 



*f 7 Q 7^ fi v * + * id v 7?^7 Q f * 



■ 95 



* r ft" it" * + * id"^ f * 



(37) 



r 



The first line in Eq. (37) is just the usual kinetic RHD con- tribution to T^, while the additional kinetic terms originate 
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from the evaluation of the higher-order derivatives in Eq. (12). 
These terms will be important for the the thermodynamic con- 
sistency of the model and the validation of the Hugenholtz- 
Van Hove theorem [14, 30]. Again the terms not shown in 
Eqs. (36) and (37) describe the contribution of terms contain- 
ing the derivatives of the meson fields. 



IV. RMF APPROACH TO INFINITE NUCLEAR MATTER 

In the mean-field approximation the mesons are treated as 
classical fields. Infinite nuclear matter is described by a static 
homogeneous, isotropic, spin and isospin-saturated system of 
protons and neutrons. In this case, the spatial components of 
the Lorentz-vector meson fields vanish with lu^ — > (cj°, 0), 
and in isospin space only the neutral component of the isovec- 
tor fields survive, i.e., p^ — > (p®, ) and S — > (S3. For sim- 
plicity, we denote in the following the third isospin compo- 
nents of the isovector fields as p and S. 

The derivation of the RMF equations starts with the usual 
plane wave ansatz 



*i(s,p) = Ui(s,p)e ip x " , 



(38) 



where i stands for protons (i = p) or neutrons (i = n) and 
p^ = (E,p) is a single nucleon 4-momentum. The applica- 
tion of the non-linear derivative operator T> to the plane wave 
ansatz of the spinor fields results in 



V(^)% = V(Oui(s,p)e-n 



(39) 



(40) 



where the regulators in the r.h.s. of above equation are now 
functions of the scalar argument £ = — Va ? . 

With the help of Eqs. (38) and (39) one gets the Dirac equa- 
tion similar to Eq. (28) with selfenergies given by 



^ n = g^V + g p T lP ^V 



S s? = g a <jV + gsnSV 



(41) 



(42) 



where now t; = +1 for protons (i = p) and n = — 1 for 
neutrons (i = n). We note again that in the RMF approx- 
imation to infinite matter the additional terms including the 
meson field derivatives vanish. This largely simplifies the for- 
malism, since these terms which show up in the original Dirac 
equation, see Eq. (28), do not appear any more. 
The solutions of the Dirac equation takes the form 



Ui(s,p) = Ni 



a ■ p 

V E* + m* ifi 



) 



with spin eigenfunctions ip s , the in-medium energy 



E* 



E 



(43) 



(44) 



and the Dirac mass 



???: 



(45) 



For a given momentum the single particle energy E is ob- 
tained from the in-medium on-shell relation 



Ef 



V 



... 2 



(46) 



The factor Ni is determined from the normalization of the 
probability distribution, that is J d 3 x J° = 1. In the con- 
ventional RMF models the baryon density is given by the fa- 
miliar expression J° = an( j the normalization condition 

/ d 3 x ft* = 1 would result in N, = ^ E " 2E T' ■ In the NLD 
model one has to use Eq. (9) for the Noether current by keep- 
ing in mind that the infinite series of meson field derivatives 
vanish in the RMF approach to nuclear matter. In this case, 
see again Appendix E for details, the conserved baryon cur- 
rent J p is resummed up to infinity and the result reads 



=E 

i—p,n 



(47) 



+g*(*i[dgV]* i )<T - g u (Vi[d£V]y a iI> 
- g fi T i (¥ i [djp]>f l * i )p a + gsntVildgV^AS 



The 0-component of the Noether current describes the con- 
served nucleon density ps — J°, from which also the relation 
between the Fermi momentum pp and ps is uniquely deter- 
mined. In particular, using the Gordon identity and Eqs. (41) 
and (42) for the RMF selfenergies, one obtains 



J" 



(27r) a 



E 



d 3 P N? 



P\<PFi 



E 



E* 



where k — 2 is a spin degeneracy factor, pp i stands for the 
proton or neutron Fermi-momentum and the effective momen- 
tum is given by 



Pi =P' 

One defines now a new in-medium 4-momentum Ylf as 



and arrives to the following expression 



/3 



E 



(2tt) 3 ^ 



d 3 P N- 



nf 

E* 



(48) 



(49) 



(50) 



\<PFi 



On the other hand, the general definition of the baryon current 
results from the covariant superposition of all the occupied 
in-medium on-shell nucleon positive energy states up to the 
proton or neutron Fermi momentum [7] 



(2^)3 ^ 

i—p.n 1 1 ^ 

\P\<PFi 



^HpTp*u 



m 



d 4 p 
2 ) 20(p°). 



(51) 



7 



In the NLD approach the mean-field selfenergies depend ex- 
plicitly on the single-particle momentum pP. Therefore, using 
the properties of the (5-function the time-like dp component 
can be integrated out explicitly. The result reads 



J" = 



(2tt)3 / 



E 



d 3 p 



p,n | _, , 



■ 



(52) 



Comparing Eq. (52) with the equation for the NLD current, 
Eq. (50), one gets the following result for the normalization 



N; 




(53) 



and the bilinear products between the in-medium spinors of 
protons and neutrons are given by 



Ui{p)u t {p) =J±, 
U t (p)j°Ui(p) =1 . 



(54) 
(55) 



Eq. (55) ensures also the proper normalization of the proba- 
bility distribution, i.e., J = 1. 

In our first work [25], where the non-linear derivative 
model has been proposed, the correction terms proportional 
to the partial derivatives of the selfenergies with respect to the 
single-particle momentum p^ in Eq. (49) were not taken into 
account. Even if their contributions are small at low densities, 
these terms will be included in the present calculations which 
attempt to consider also the high density domain of the EoS in 
neutron stars. On the other hand, the inclusion of these terms 
is crucial for a fully thermodynamically consistent formalism 
and is independent of the particular form of the cut-off func- 
tions. Note that, the additional cut-off dependent terms in the 
baryon and energy densities of Ref. [25] are now canceled by 
the proper normalization conditions. 

The energy-momentum tensor in NLD is obtained by ap- 
plying the Noether theorem for translational invariance. In 
nuclear matter the resummation procedure results in the fol- 
lowing expression 



-gT{C). 

The evaluation of the expectation values in Eq. (56) can be done in a similar way as for the current with the result 



\p\<PF i 



* Li , * 

Pi +rn i 



-gT(C) 



(56) 



(57) 



Using Eq. (49) one arrives to the final expression for the 
energy-momentum tensor in the NLD formalism, which can 
be written in the following form 



P (2tt)3 J p n: 



■^-g^(C), (58) 



from which the energy density e = T 00 and the pressure P 
can be calculated, i.e., 



E 



(2tt) 3 



d 3 pE(p)-{C) 



|p|<p* 



(2tt) 3 

\P*\<PF i 

Eqs. (59) and (60) look similar as the familiar expressions of 
the usual RMF models. However, the non-linear effects in- 



duced by the regulators show up through the generalized mo- 
mentum IT' 1 and through the dispersion relation for the single- 
particle energy E(p). Note the different form of the gener- 
alized momentum, II' 1 , when one chooses energy or momen- 
tum dependent cut-off functions. Indeed, in the latter case the 
spatial derivatives in IT contribute in the pressure, while for 
energy-dependent cut-off functions they vanish and 11, = p 
holds. In any case, the expressions for the energy-density and 
pressure within the conventional RMF models are recovered 
by simple replacement IlJ 1 — > p*^, which is just equivalent to 



(59) taking the limiting case A — > oo in the NLD expressions. 



Finally, the NLD meson-field equations in the RMF ap- 
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proach to nuclear matter read 

dU 



9cPs , 



(61) 



m 2 ,uj 



(*i7°^*i) = 9uPo , (62) 



i— p,n 



™ P P =.9 P 



^ T-i/tfiT ©**) = fl p pj , (63) 



m]5 =g s J2 rJ^iV^i) = 



gspis 



(64) 



i— p,n 



Using Eqs. (54) and (55), the evaluation of the source terms 
of the meson-field equations is straightforward. In particular, 
the scalar-isoscalar p s , vector-isoscalar pg, vector-isovector 
pi and scalar-isovector pj$ are given 

P' = i^ E f**m v v>> (65) 

|p|<p F . 



PO 



(2tt) 3 ^ 

|p|<pj?. 



P/ — POp - POn ) 



PIS = Psp — Psn ■ 



d 3 p^V(p), 



(66) 



(67) 



(68) 



The meson-field equations of motion show a similar structure 
as those of the standard RMF approximation. For example, 
the scalar-isoscalar density p s is suppressed with respect to 
the vector density po by the factor m* /II?, in a similar way 
as in the conventional Walecka models [3]. However, the sub- 
stantial difference between NLD and conventional RMF ap- 
pears in the source terms which now contain in addition the 
momentum-dependent regulator V. 



V. RESULTS 

A. Model parameters 

The non-linear derivative operators, T>, in the NLD La- 
grangian are not constrained from first principles and allows 
us to consider different functional forms of T>. In nuclear mat- 
ter, these regulators can be chosen as functions of the single- 
particle energy or momentum, depending on the choice of 
the auxiliary multi-dimensional parameters £f . The available 
constraints are the bulk properties of nuclear matter and the 
empirically known energy dependence of the in-medium opti- 
cal potential. It is also well known [9, 31, 32] that the selfen- 
ergies should decrease or saturate as function of baryon den- 
sity and single-particle 4-momentum. In the NLD model all 
these features of the relativistic mean-fields can be realized us- 
ing energy or momentum dependent form factors which reg- 
ulate the high energy (momentum) behavior of the nucleon 
4-momentum. 



Various choices of the regulator functions are possible. We 
have done calculations for different forms of V and found 
that the simplest momentum dependent monopole form fac- 
tor provides the best description of the low and high density 
nuclear matter properties and agrees very well with the empir- 
ical momentum dependence of the in-medium Schrodinger- 
equivalent optical potential. A momentum-dependent cut-off 
of a monopole form can be obtained by the following choice 



I] 7 = 1 



(69) 



with u? = (0, 0, 0, 0), v% = (0, 1, 0, 0), = (0, 0, 1, 0) and 



(0, 0, 0, 1). In nuclear matter this results in 



V 



A 2 



A 2 



■P' 



(70) 



Furthermore, in our fit to bulk properties of nuclear matter we 
use different cut-off parameters A s = A a and A„ = A u = A p 
for the scalar and vector meson-nucleon vertices, respectively, 
and we neglect in the following the contribution of the 5- 
meson. 

The NLD model contains in total eight parameters. These 
are the meson-nucleon couplings g a , g u and g p , the parame- 
ters b and c of the selfinteractions of the er-meson, the mass 
m a of the cr-meson, and the cut-offs A s and A„. In principle, 
the masses of the uj- and p-mesons should be also included 
in the fit. In all the calculations concerning the fit the results 
for these two masses turned out to be always around their free 
values. Therefore, we keep for m u and m p the bare masses. 
Seven parameters, i.e., g a , g u , b and c, A s , A.„ and m a are 
adjusted to the bulk properties of symmetric nuclear matter 
and to the empirical energy dependence of the in-medium op- 
tical potential. The remaining parameter g p is determined by 
the experimentally known symmetry energy value at satura- 
tion density. The constraints in symmetric nuclear matter are 
the binding energy per nucleon Eb = e / ps — m and the com- 
pression modulus K = Qp 2 sat % ~f b \p B =p 3at at the ground-state 
or saturation density p sa t, and the saturation density p sat it- 
self. Furthermore, the momentum dependence is fixed by the 
optical potential U op t at ground-state density and at two ki- 
netic energies. First, at Eki n = 205 MeV (where U op t = 
MeV) and at E km = 1000 MeV. 

The numerical calculations for the fit procedure have been 
performed using the Nelder-Mead minimization algorithm 
NELMIN [33]. Furthermore, an adaptive non-linear least- 
square package NL2SOL [34, 35] has been supplemented in 
order to test robustness of the fit procedure. The experimen- 
tal data used in the fit are supplemented by the correspond- 
ing errors as provided by the empirical information in Table I. 
At each iteration step of the minimization routines, the cou- 
pled set of the NLD equations has been solved with a Powell- 
hybrid method as provided by the HYBRD routine [36]. The 
integrals which appear in the source terms of the meson- 
field equations have been treated numerically using an adap- 
tive Gauss algorithm. Note that, in the NLD model with the 
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Model 


psat 


E b 


K 




L 






Ref. 




[fm- 3 ] 


[MeV/A] 


[MeV] 


[MeV] 


[MeV] 


[MeV] 


[MeV] 




NLD 


0.156 


-15.30 


251 


30 


81 


-28 


-514 


this work 


NL3* 


0.150 


-16.31 


258 


38.68 


125.7 


104.08 


-650.12 


[37] 


DD 


0.149 


-16.02 


240 


31.60 


56 


-95.30 


-431.30 


[19] 


D 3 C 


0.151 


-15.98 


232.5 


31.90 


59.30 


-74.7 


-430.50 


[22] 


DBHF 


0.185 


-15.60 


290 


33.35 


71.10 


-27.1 


-453.70 


[38, 39] 




0.181 


-16.15 


230 


34.20 


71 


87.36 


-340 


[40] 


empirical 


0.167 ±0.019 
[41^3] 


-16 ± 1 
[41-43] 


230 ± 10 
[43^15] 


31.1 ± 1.9 
[46] 


88 ± 25 
[46] 




-550 ± 100 
[47] 





TABLE I. Bulk saturation properties of nuclear matter, i.e., saturation density p sa t, binding energy per nucleon Eh, compression modulus 
K, asymmetry parameter a sym , slope and curvature parameters L and K 3ym , respectively, and the observable K a3y in the NLD model. Our 
results are compared with the non-linear Walecka parametrization NL3*, the density dependent DD and the derivative coupling D 3 C models 
as well as with two versions of the microscopic DBHF approach. The empirical values are shown too. 







cut-off 


A s A v 
[GeV] [GcV] 


9o g P b c 
[fm- 1 ] 


m CT m„ m p 
[GeV] [GeV] [GeV] 


NLD 


1 


A 2 
A 2 +p 2 


0.95 1.125 


10.08 10.13 3.50 15.341 -14.735 


0.592 0.782 0.763 



TABLE II. The parameters of the NLD model. First and second columns show the form of the non-linear operator and the regulator in nuclear 
matter, respectively. The other columns show the values of the parameters, as resulted from the fit to nuclear matter bulk properties. 



momentum-dependent regulator functions the solution of the 
dispersion relation for the single particle energy E, Eq. (46), 
does not involve additional complex root-finding algorithms, 
because the RMF selfenergies depend on the nucleon momen- 
tum p only, and not on the particle energy as in our previous 
works [23-25]. 

The functional form of the non-linear operator, its regula- 
tor (cut-off) in nuclear matter and the fit parameters of the 
NLD model are shown in Table II. Table I shows the ex- 
tracted bulk properties of nuclear matter, i.e., the binding 



energy, the compression modulus, the asymmetry parameter 
a sym = E sym (p sat ) (E sym is the symmetry energy), the 
slope and curvature parameters, L = Spsat ^^ \ PB=p , at 

andK asy = K sym - QL (with K sym = 9p 2 sat 9 g^ m \ PB=Psat 
being the compressibility of the symmetry energy), respec- 
tively, in the NLD model, and in comparison with other RMF 
models widely used in the literature. These are the NL3* 
parametrization [37], the density dependent (DD) and deriva- 
tive coupling (D 3 C) models [19, 22]. The bulk nuclear matter 
properties in the NLD model are comparable with those re- 
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P B [fm" 3 ] 

FIG. 1 . Vector (upper panel) and scalar (lower panel) selfenergies as 
function of the baryon density ps in the NLD model (solid curves) 
and in the DBHF approach (filled squares) [39]. The calculations 
refer to isospin-symmetric (a = 0) nuclear matter. 

suits in the NL3*, DD and D 3 C parameterizations, while for 
the saturation density and the slope and curvature parameters 
of the symmetry energy the NLD calculation is closer to the 
empirical data. It is interesting that the NLD results are also 
comparable and in some cases even better than the DBHF cal- 
culations. 

As we will discuss later on in more details, the NLD model 
with the same parameters of Table II describes remarkably 
well the empirical energy dependence of the Schrodinger- 
equivalent optical potential. This is not the case in standard 
RMF approaches, such as the NL3* and DD models, except 
if one uses supplementary derivative interactions (D 3 C), how- 
ever, with the cost of many additional parameters [22]. 



B. The NLD equation of state 

We start the discussion of the NLD calculations with the 
density dependence of the relativistic mean-fields. Fig. 1 
shows the density dependence of the nucleon selfenergies for 
isospin-symmetric nuclear matter in the NLD model. Our re- 
sults are compared also with DBHF calculations, widely used 
in the literature [39]. It is remarkable that a simple folding of 
the meson-nucleon couplings gi (i = <t,uj, p) with the cut-off 
function V results in a very similar density dependence of the 
of NLD selfenergies as compared to the DBHF selfenergies, 
in particular, for densities above pb > 0.3 — 0.4 /m~ 3 . Note 
again that for momentum dependent cut-off functions the stan- 
dard normalization of the spinors, i.e., Ni = y5JpS is not 

affected in this case, since II? = E* holds, and according to 
Eq. (49) only the pressure is modified, as shown in Eq. (60). 
Fig. 2 shows the equation of state in terms of the binding 




-50 1 ■ 1 ■ 1 ■ 1 ■ 1 ■ ' — 

0.2 0.4 0.6 0.8 1 

P B [ fm " 3 l 

FIG. 2. Equation of state (EoS) in terms of the binding energy per 
nucleon as function of the baryon density ps for isospin-symmetric 
(a = 0, lower curve and symbols) and pure neutron matter (a = — 1 , 
upper curve and upper symbols). The NLD results (solid curves) are 
compared with DBHF calculations (filled squares) [38]. 



energy per nucleon as function of baryon density for isospin- 
symmetric (a = 0) and pure neutron matter (a = —1). The 
isospin-asymmetry parameter a is defined in the usual way 

J°—J° 

as a — j p u +J Z , where J°„ denote the proton and neutron 
densities, respectively. The momentum-dependent monopole 
form-factor of the NLD model regulates the high-density de- 
pendence of the fields such that the NLD EoS fits very well 
the DBHF calculations, for both, symmetric nuclear and pure 
neutron matters. Note that the NLD parameters are not fit- 
ted to the calculations of DBHF models, but to the empirical 
information at ground state density only. 

In Fig. 3 the density dependence of the EoS is displayed 
again, but now we compare the NLD calculation with other 
RMF models, which have been widely used in studies of finite 
nuclei and of nuclear matter. Even if all models give similar 
results for the saturation point of nuclear matter, the differ- 
ences between the NLD model and the other RMF approaches 
at densities beyond the ground state are large. In particular, 
the NL3* parametrization of Lalazissis et. al. [37], which is 
based on the Walecka model with non-linear selfinteractions 
predicts the stiffest EoS, while the RMF approaches DD and 
D 3 C with density-dependent meson-nucleon couplings give a 
softer density behavior at high densities. However, none of 
these RMF approaches can reproduce the microscopic calcu- 
lations of the DBHF model at such high densities. In gen- 
eral, the non-linear density dependence of the NLD model, 
induced by the cut-off form factor, result in a much softer EoS 
at high densities, and it is in agreement with the DBHF calcu- 
lations. Note that heavy-ion studies at energies just below the 
kaon production threshold (Ei, eam < 1 GeV/A) on collec- 
tive nucleon flows [48-50] and on produced mesons (positive 
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FIG. 3. Same as in Fig. 2, but now the NLD model (solid 
curve) is compared with various RMF calculations. These are 
the NL3* parametrization (dashed-dotted-dotted curve) [37], the 
density-dependent (DD, dashed curve) and derivative-coupling 
(D 3 C, dashed-dotted curve) models [19, 22] and the linear Walecka 
model (RHD, dotted curve) [3]. The DBHF calculations (filled 
squares) are shown too. All the calculations refer to isospin- 
symmetric (a = 0) nuclear matter. 



charged kaons) [51-54] support a rather soft EoS at densities 
around p B ~ (2 - 2>)p sat . 

Furthermore, an interesting issue concerns the inclusion of 
the non-linear selfinteractions of the cr-field in the RMF de- 
scriptions. It is well known, that such terms may cause diver- 
gences at very high densities, where one intends to study com- 
pact neutron stars. Indeed, many RMF approaches such as 
the NL3* [37] or the NLp and NLpS [55] parameterizations, 
lead to a non-physical behavior of the cr-field at supra-normal 
densities. This is not the case in the NLD model, where the 
non-linear selfinteraction terms of the er-meson does not cause 
any divergences of the cr-field even at supra-normal densities. 
This novel NLD feature arises from the suppression of the 
scalar sources by the cut-off function T>, which in the limiting 
case of large densities ps tend to vanishes. 

Another quantity of interest is the symmetry energy. This 
quantity is important for astrophysical applications, since 
it directly influences the density dependence of the proton- 
fraction in the neutron star interior. It is extracted as the dif- 
ference between the EoS's from pure neutron and symmetric 
matters. Fig. 4 shows the symmetry energy as function of 
the baryon density in the NLD model, in other RMF mod- 
els and in the DBHF approach. Note that the NLD results 
differ from the previous works [23]. This is because now 
the NLD model consistently includes the renormalization of 
the spinors. This is essential for the thermodynamic consis- 
tency of the formalism, which was not taken into account in 
Ref. [23]. Furthermore, a different form of the regulator is 
used here and the parameters of the NLD model are fitted us- 



FIG. 4. Same as in Fig. 3, but now for the symmetry energy E sym as 
function of baryon density ps- The additional symbols (diamonds, 
triangle and circles) around and below the saturation density refer to 
empirical data [56, 57], 

ing low-density observables at saturation. At first, all models 
describes fairly well the empirically known region around the 
saturation density, as it was already shown in Table I. How- 
ever, the differences between the models become more pro- 
nounced with increasing density. The most stiffer symmetry 
energy is obtained again in the NL3* and the RHD calcula- 
tions, because in these cases E sym is just proportional to the 
p-meson, which linearly increases with density. Furthermore, 
the symmetry energy is considerably softened in the DD and 
D 3 C approaches, due to the exponentially decreasing density 
dependence of the p-nucleon vertices. On the other hand, the 
NLD model leads to a softer density behavior of the symme- 
try energy relative to the standard NL and RHD parametriza- 
tion but to a stiffer dependence than the density-dependent 
approaches. It is again an interesting feature that the NLD 
calculations predict almost perfectly the calculations of the 
microscopic DBHF calculations [38]. 

C. In-medium nucleon optical potentials 

Apart the density dependence, which arises from both, the 
source terms of the meson-field equations and the cut-off 
functions, the NLD selfenergies depend explicitly on momen- 
tum. The momentum (or energy) behavior of the RMF fields 
is described by the in-medium Schrodinger-equivalent optical 
potential. This quantity results from the reduction of the Dirac 
equation and reads [58] 

U opt = -S+-V+^-(S 2 -V 2 ) , (71) 

where the selfenergies S = ^ S i(pB,p), V = ^vi(Pb,p) re- 
fer to a proton (i = p) or a neutron (i = n) with a particular 
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FIG. 5. (Color online) The in-medium proton optical potential U op t 
according Eq. (71) as function of the in-medium single-particle ki- 
netic energy, Eq. (72), within the NLD model (solid curve) and the 
the RHD [3] approach at saturation density. The symbols (filled cir- 
cles) refers to the results of the Dirac phenomenology [5, 6]. 




FIG. 6. Energy dependence of the Lane-type optical potential Ui so , 
Eq. (73), for asymmetric (a = 0.4) nuclear matter at saturation 
density. Calculation in the RHD [3] (dotted curve) and NLD (solid 
curve) models are shown. 



momentum \p \ = p relative to nuclear matter at rest at a given 
density ps and isospin-asymmetry a. The single-particle en- 
ergy E is calculated from the in-medium dispersion relation 

E = Vm* 2 + p 2 + V (72) 

and the kinetic energy reads Ekm = E — m [9, 22, 31]. Here 
we consider the real part of U op t only, which enters the mo- 
mentum dependent mean-field dynamics of a nucleon in nu- 
clear matter. The imaginary part, which is beyond the scope 
of the present work, can be further accounted for in the colli- 
sion integral within the relativistic transport equation, see for 
instance Ref. [12], when applying the NLD approach to the 
description of heavy-ion collisions. Alternatively, one can use 
a dispersion theoretical framework to calculate the imaginary 
part of the optical potential using as an input the real part of 
the optical potential in the dispersion integral, see Ref. [24] 
as an example of such approach to the imaginary part of U opt- 
In this context we would like to note, that the imaginary con- 
tributions of the relativistic mean-fields do not influence es- 
sentially the real part of the optical potential, see for exam- 
ple Refs. [7, 22]. The reason is, that the imaginary contribu- 
tions of the selfenergies, as obtained from empirical studies, 
are of similar magnitude. They enter only through the terms 
quadratic in the selfenergies in the expression for the real part 
of the Schrodinger-equivalent optical potential, Eq. (71), and 
therefore they almost cancel each other. 

Fig. 5 shows the results of our calculations in comparison 
with the empirical data extracted from Dirac phenomenol- 
ogy [5, 6] for the in-medium proton optical potential in sym- 
metric (a = 0) nuclear matter at saturation density pb = Psa t - 
The optical potential saturates with increasing single-particle 
energy in the NLD model and it is in agreement with the 



experimental data. The saturation mechanism is induced by 
the explicit momentum dependence of the cut-off functions, 
which drop with increasing momentum. 

It is a novel feature of the NLD model, that the regula- 
tors with a simple momentum-dependent monopole form are 
sufficient to describe accurately the bulk properties of nu- 
clear matter and at the same time the empirical energy de- 
pendence of the optical potential. In fact, this issue has been 
a long-standing problem in nuclear matter studies, when one 
attempted to describe heavy-ion reactions within RMF mod- 
els [59]. As also shown in Fig. 5, in standard RMF models, 
such as the widely used NL-parameterizations [37, 60], which 
describe excellent the saturation properties and also the prop- 
erties of finite nuclei, cannot reproduce the correct energy de- 
pendence of the optical potential. Moreover, they strongly 
diverge with increasing energy. Similar conclusions are ob- 
tained in the RMF approaches with density-dependent ver- 
tex functions [19], except if one includes additional energy- 
dependent terms [22]. Note that, the microscopic DBHF mod- 
els [9, 31, 32] describe satisfactory the empirical data at low 
energies bellow the pion production threshold. 

For isospin-asymmetric nuclear matter the key quantity is 
the so-called Lane potential [46, 61], which is defined as the 
difference between the neutron and proton optical potentials 

U lso = " p . (73) 
2\a\ 

Fig. 6 displays the energy dependence of the Lane potential. 
In contrast to the case of isospin-symmetric matter, empirical 
information is here less known. The studies from Ref. [46] 
predict a decrease of the Lane potential with increasing mo- 
mentum, but other analysis [61] conclude the opposite trend 



13 




-100 



250 



500 750 
E„„ [MeV] 



1000 



1250 



FIG. 7. Energy (momentum) dependence of the optical U op t (Lane) 
potentials in the main (inserted) panel. The curves have the same 
meaning as in Figs. 3 and 4. Here the D 3 C results are taken 
from [64]. 



with the result of an increasing Lane potential. The micro- 
scopic DBHF calculations predict also a decreasing optical 
potential, which is in agreement with results of other BHF 
calculations [62, 63]. Furthermore, the standard RHD model 
leads to an almost quadratic (or linear in energy) dependence 
in momentum, just because of the missing momentum depen- 
dence in the RHD selfenergies (similar results are obtained in 
the NL* and DD parameterizations). The NLD calculations 
predict a weakly decreasing Lane potential, which is in quali- 
tative agreement with the (D)BHF results. 

We compare now the NLD results separately in Fig. 7 for 
both, the in-medium proton optical potential for symmetric 
nuclear matter and the Lane potential, with the same RMF 
approaches as in Figs. 3 and 4. Not only the original lin- 
ear Walecka model (RHD), but also more modern RMF ap- 
proaches, such as the non-linear NL3 model in its updated 
version (NL3*) and the density-dependent RMF model (DD) 
predict an energy dependence, which is not consistent with 
phenomenology (filled circles). This circumstance was im- 
proved in a modified DD model (D 3 C) by the inclusion of 
terms linearly proportional to the single-particle energy. The 
NLD model provides here a smoother energy dependence. 

The situation for the Lane potential is shown in the insert 
of Fig. 7. The interesting issue here is, that not only the stan- 
dard RMF models (RHD, NL3* and DD), but also the energy- 
dependent RMF approach (D 3 C), show a common behavior in 
momentum. This is due to the fact that the isovector channel 
in the D 3 C RMF model does not include any explicit momen- 
tum dependence [64]. On the other hand, in the NLD model 
also the isovector mean-field depends on momentum and pre- 
dicts a decrease of the Lane potential with increasing momen- 
tum. This NLD trend seems to be supported by microscopic 
DBHF (filled symbols in the insert of Fig. 7), at least on a 
qualitative level. 
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FIG. 8. Thermodynamic consistency of the NLD model for pure 
symmetric nuclear (a — 0, left panel) and asymmetric (a = —0.4, 
right panel) matters. The pressure densities P as function of the 
baryon density pa are shown, as calculated within the perfect fluid 
formula (solid curves) using the l.h.s. of Eq. (76) and from the ther- 
modynamic definition (filled diamonds) using the r.h.s. of Eq. (76). 



D. Thermodynamic consistency of NLD model 

As showed by Weisskopf [30] in an independent-particle 
model, the single particle energy at the Fermi surface must 
be equal the average energy per particle at saturation density. 
Hugenholtz and Van Hove [14] proved this also for an inter- 
acting Fermi gas at zero temperature. 

For the thermodynamic consistency of the RMF model it is 
sufficient to prove Euler's theorem 



e= -P 



E 

i— p : n 



(74) 



with the chemical potentials /Ltj = 
namic definition of the pressure 

2 d(e/ PB ) 
" — Pb~ 



dpi 



and the thermody- 



(75) 



The expression in Eq. (74) allow to examine the internal 
consistency of the model. For this purpose it is sufficient to 
check the equality between the pressure obtained from the as- 
sumption of nuclear matter as a prefect-fluid system and the 
pressure obtained from the thermodynamic definition. That is 



P 



/rjixx _|_ rpyy -|- J<ZZ\ 



2 d(e/p B ) 
Pb~ 



dpB 



(76) 



where e and T" (i = x, y, z) denote the energy density and the 
spatial diagonal components of the energy-momentum tensor 
tensor T^ v , respectively. 
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FIG. 9. Equation of state in terms of the pressure densities as function 
of the baryon density (normalized to the corresponding saturation 
values p S at)- The shaded areas denote possible experimental regions, 
as extracted from studies on heavy-ion collisions [48]. The different 
curves have the same meaning as in Figs. 3 and 4. 



Conventional RMF approaches with bare meson-nucleon 
couplings or density-dependent meson-nucleon vertices are 
thermodynamically consistent [4, 21, 60, 65]. However, in the 
case of explicit energy or momentum-dependent mean-fields 
the situation is more complex. In fact, as we have examined 
here, one has to take care of the proper renormalization of the 
Dirac fields. This issue was not taken into account in previ- 
ous studies [25] resulting in a deviation between the l.h.s. and 
r.h.s. of Eq. (76) by ~ 5 — 10% at high baryon densities. 

We have checked that the NLD model presented here sat- 
isfies thermodynamic consistency exactly for any kind of 
energy- or momentum-dependent cut-off form factors. This is 
demonstrated in Fig. 8 for the monopole form factor, where 
the pressure as function of density for symmetric nuclear 
matter (left panel) and pure neutron matter (right panel) are 
shown. Indeed, the pressures obtained from both definitions 
of Eq. (76) agree exactly. Furthermore, the chemical poten- 
tials at zero temperature are equal to the corresponding Fermi 
energies of protons and neutrons (not shown here) and there- 
fore Euler's theorem is evidently fulfilled. 

We can now consider the high-density behavior of the pres- 
sure by comparing our results with available empirical infor- 
mations. Fig. 9 shows the density dependence of the pressure 
densities in symmetric matter and pure neutron matter. The 
conventional RHD approach leads again to a stiff behavior for 
all the densities. Similar conclusions are drawn for the NL3*, 
DD and D 3 C approaches for symmetric matter, where for pure 
neutron matter the density-dependent models come closer to 
the empirical HIC data. A more detailed discussion on these 
approaches can be found in Ref. [66]. The pressure in the 
NLD model show generally a softer density dependence and 
agree better with the estimated experimental regions, in par- 



FIG. 10. Relation between Neutron star mass M (in units of the so- 
lar mass Mq) versus radius R in the NLD model (dashed and solid 
curves). The three horizontal shaded bands refer to astrophysical 
measurements from double neutron star (NS) systems [73, 74] and 
from the pulsars PSR J1903+0327 [75] and PSR J1614-2230 [67], as 
indicated. The other shaded areas bordered by thick curves indicate 
parameter space excluded by general relativity, causality (shaded 
area on the top-left) and rotational constraints (shaded area on the 
bottom-right) [73, 74]. 

ticular, at densities up to ps — ^Psat for the symmetric case 
and for all densities for pure neutron matter. Note that for 
larger densities conclusions on the nuclear matter EoS from 
heavy-ion studies are more ambiguous, because at such high 
densities (or corresponding beam energies larger than 4 GeV 
per nucleon) a large fraction of the initial energy is converted 
into new degrees of freedom [48, 66]. 



E. High density observables and Neutron Stars 

We test now the high density domain of the NLD equa- 
tion of state. Compact neutron stars offer such an oppor- 
tunity to gain deeper insight into compressed baryonic mat- 
ter. Of particular interest are recent measurements on the bi- 
nary millisecond pulsar J1614-2230 with a mass of (1.97 ± 
0.04) Mq [67]. The latter is much heavier than the average 
mass of the binary radio pulsars M = 1.35 ± 0.04 Mq [68] 
and provides a strong constraint on high density EoS. There- 
fore, we apply the NLD model to spherical, non-rotating 
stars in /3-equilibrium between protons, neutrons and elec- 
trons including crustal effects on the EoS [69]. The star struc- 
ture is calculated by solving the Tolman-Oppenheimer-Volkov 
(TOV) equation [70-72]. 

The results for neutron stars are shown in Fig. 10 in terms 
of the mass-radius relation. The various astrophysical mea- 
surements of NS masses [73, 74] can be arranged in the three 
horizontal shaded areas as displayed in Fig. 10. The low- 
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FIG. 11. Pressure densities as function of the baryon density for nu- 
clear matter in /^-equilibrium within the NLD (solid curve) model. 
The shaded areas together with the three error bars indicate the 
pressure-density region extracted directly from neutron- star mea- 
surements [26]. 



est band around an average value of 1.44 Mq refer to the 
well established measurements on double neutron star sys- 
tems and the middle one around 1.67 M Q on the extracted 
mass of the pulsar PSR J1903+0327. The band on the top 
represents by far the highest precisely observed neutron star 
mass 1.97 ± 0.04 M of the pulsar PSR J1614-2230 [67]. 
There are two regions in Fig. 10 excluded [73, 74] by general 
relativity, causality (shaded area on the top-left) and rotational 
constraints (shaded area on the bottom-right). 

The neutron stars mass-radius relation in the NLD model is 
shown by the solid/dotted curve in Fig. 10. The dotted part 
of the NLD curve is excluded by rotational constraints. The 
solid curve crosses the low-mass constraints, and arrives to a 
maximum neutron star mass of M = 2.03 M at a radius of 
R = 11.07 km and a corresponding central density of p c ~ 
7 p S at- The NLD prediction for the maximum value of neutron 
star masses crosses also the constraint provided by the pulsar 
PSR J1614-2230, and therefore this recent mass measurement 
is accommodated by the NLD model. 

Other possible constraints on the high-density EoS are ob- 
tained by statistical Bayesian analyses, which rely on neutron 
star measurements [26, 27]. They provide the most proba- 
ble distribution of the equation of state, as shown in Fig. 1 1 
(shaded area) for highly compressed matter in /^-equilibrium. 
The RHD model (not shown in this figure) leads to a too stiff 
density dependence and it overestimates this empirical region 
particularly at high densities. The NLD calculations, where 
only nucleonic degrees of freedom are accounted for, describe 
fairly well the most probable region of the pressure at high 
densities. Note again that the parameters of the NLD model 
have been adjusted just to the saturation properties of nuclear 
matter and to the energy dependence of the in-medium pro- 



ton optical potential at saturation, without the consideration 
of any other high density observables in the fit procedure. We 
conclude here that the NLD describes well the available low- 
and also high-density constraints on EoS of nuclear matter and 
neutron stars. 



VI. SUMMARY 

In summary, in the present work the generalized form of 
the energy-momentum tensor in the NLD model was derived 
which allowed us to consider different forms of the regulator 
functions in the NLD Lagrangian. The thermodynamic con- 
sistency of the NLD model was further demonstrated for ar- 
bitrary choice of the regulator functions. A thorough study of 
the properties of nuclear matter around saturation density has 
been performed. We have shown that the NLD approach de- 
scribes well the saturation properties of the nuclear matter and 
compares remarkably well with microscopic calculations and 
Dirac phenomenology. We have investigated the high density 
part of the NLD EoS. This is relevant for the neutron stars in 
/3-equilibrium. We found that the low density constraints im- 
posed on the nuclear matter EoS and by the momentum depen- 
dence of the Schrodinger-equivalent optical potential lead to a 
maximum mass of the neutron stars around M ~ 2M . The 
latter mass accommodates the observed mass in the J1614- 
2230 millisecond radio pulsar system. We further studied the 
EoS of matter in /3-equilibrium and find that the high den- 
sity pressure-density diagram as extracted from astrophysical 
measurements can be well described in the NLD model which 
rely on nucleonic degrees of freedom only. 

The EoS proposed here can be used in transport theoretical 
studies of nuclear collisions, since, it describes very well both, 
the low energy (density) and the high energy (density) regions 
of the nuclear phase diagram. The model predicts saturation 
of the optical potential at high energies and results in satura- 
tions of the symmetry energy. An interesting finding is that 
the momentum dependent interaction make the EoS softer at 
low densities, however, it is still stiff enough at supra-normal 
densities to account for the recent measurements of the neu- 
tron star masses. Furthermore, the model can be applied to 
the transport description of the anti-nucleon optical potential 
as well as to the study of dynamics of compressed matter in 
reactions induced by heavy-ions and anti-proton beams at the 
future FAIR facility. 
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Appendix A: Notations and inline tisemal variations 

We use following abbreviations for higher-order partial 
derivatives 

ai ■■■a n ■ — 9 ai a2 ■ ' • 9a n 3 
QOci—a n ■—Q a iQ a 2...Q a n 

In the following appendices we will need various defini- 
tions of infinitesimal variations, which are specified here. The 
total variation of a field <p r (x) with x M = (x°,x ) is defined 

as 

8 T <p r {x) := v' r {x') - ip r (x) , (Al) 
with the variation with respect to the 4-coordinate given by 

x'» := x» + 6x» =x^ + M» v x v + e» , (A2) 

where Sx^ denotes an infinitesimal transformation, e.g., a con- 
stant translation, and/or a rotation, A-d^x^, fi^ is an in- 
finitesimal antisymmetric tensor and A??'"' = — A$ I//i . We 
define the infinitesimal transformation of the field at a fixed 



argument as 

Sip,- = (fi' r - (f r • (A3) 

For the derivation of the Noether theorem we will need not 
only the infinitesimal variation of a field at fixed argument 
only, but also of its higher-order derivatives, i.e., S(d ai (p r ), 
5(d ai a 2 <p r ), • • ■ , 5(d ai ---a n ( fr)- For such variations the com- 
mutation between symbols <5 and d is holds, that is 

5 (d ai tp r ) =d ai ip' r (x) - d ai ip r (x) 

=8 ai (lfi' r ( x ) ~ 'Prix)) = d ai Slfi r (x) , 
S (d aia2 ip r ) =d aia2 ip' r (x) - d aia2 lfi r (x) 

=<9 QlQ2 (<p' r (x) - ip(x)) = d aia2 (Sip r (x)) , 

and obviously for higher-order fields. The total variation of a 
field ifir(x) can be written in a more handleable way as 

8 T ip r (x) :=ip' r (x') - ip r (x) 

= [<p' r (x') - <p r (x')] + [<p r (x') - <p r (x)] , (A4) 

where the first term is just the variation 6(p r at fixed argument. 
The second term in Eq. (A4) is the variation with respect to the 
argument. Eq. (A4) reduces in first order to 

Sr^rix) — 5ip r (x) + d a tp r {x)8x a . (A5) 



I 

Appendix B: Derivation of the Noether theorem for higher-order Lagrangians 

For the derivation of the Noether theorem we start with the following Lagrangian density 

C = C[lfi r (x),d ai ipr(x), - ■ ■ ,d ai ... an ip r (x)] . (Bl) 

Invariance of the Lagrangian density, (Bl), with respect to an infinitesimal transformation of all the fields and their coordinates 
implies 

£[ip' r (x' ),d' ai ip' r (x'), - ■ ■ ,d' ai ... an ip' r (x')] = C [<p r (x),d ai ip r (x), ■ ■ ■ ,d ai ... an <p r (x)] . (B2) 
In terms of the total variation Eq. (B2) results in 

8 T C = £[cp' r (x'),d' ai (p' r (x'),--- ,d' ai ... otn tp' r (x')] -£[tp r (x),d ai tp r (x),--- ,d ai ... an tp r (x)} =0, (B3) 
d d 

where d u = — — and d' = , . Our goal is to arrive from Eq. (B3) to a continuity equation of the form 

ox^ F ox v 

d^F = , (B4) 

with P L being a conserved current to be determined in the following. We will work out here all the analytical evaluations up to 
third order in the partial derivatives of the Lagrangian density, and we will give all the terms up to infinity in the final equations. 
At first, the total variation of the Lagrangian density, Eq. (B3), can be written as 

S T C =C[ip' r {x'), d a tp' r (x'), d a p(p' r (x'), ...}- C[(p r (x'), d a (p r (x'), d a ptp r {x'), . . .] (B5) 
+£[ip r (x'),d a ip r (x'),d al3 ip r (x'), . . .] - C[<p r (x),d a <p r (x),d a p<p r (x), . ..] = 0, 

where we used d'ip'(x' ) = d v ip'{x' )d'x v = d v <p'(x' ) S " = d^ip'(x' ) where is a Kronecker symbol. The first line in Eq. 
(B5) is just the variation of £ with respect to the fields ip r (r'), d a cp r (r'), • ■ • at fixed argument, whereas the terms in the second 
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line give the variation of the Lagrangian with respect to the argument x^. Since we consider infinitesimal transformations only, 
it is sufficient to evaluate latter quantity up to first order with respect to the argument. Therefore, Eq. (B5) can be written as 



S T C = SC + d a CSx a = . 
The variation at fixed argument, 8C, can be evaluated in the usual way according 



dipr o{d ai <p r ) 



dC 



-S(d aia2 ip r ) + 



dC 



fr) + 



d(d aia2 if r ) - d(d a fr) 

Replacing the first term in Eq. (B7) with the help of the Euler-Lagrange equations of motion, Eq. (7), results in 



(B6) 



(B7) 



d ai 



dC 



d(d ai <p r 
dC 



dC 



d(d ai tp r ] 



6(d ai ip r ) + 



d(d aia2 tp r 
dC 



dC 



d(d aia2 (p r ) 



8{d ai a 2 fr) + 



d(d 

a i a 2 a 3 

ip r ) 

dC 



Sip r 



(B8) 



d(d tp r ) 



fr) + 



We use the commutation between the variation at fixed argument and the partial derivative (see appendix A) and apply once 
the product rule as follows for the term proportional to d aia2 5ip 



dC 



dC 



d{d 

QlQ2 fr) 1 2 1 d(d a \a 2 fr) 

For the term proportional to d aia2Ct3 6(p we apply the product rule twice 
dC dC 



d a2 5ip r + d ai 



dC 



d{d aia2 fr 



-d a2 8ip r 



(B9) 



d ai 



d a . 2 



d(d a 

\OL2OL3 fr 



d a3 5fr 

da 3 8f r 



d ai 



dC 



d{d fr) 



da 2 a 3 $fr 



(BIO) 



and so forth for the terms proportional to higher-order derivatives. In total, this procedure leads to a series of terms proportional 
to a 4-divergences only. It is more convenient to arrange these terms such to obtain several infinite series for each derivative 
field Sip, d ai Sip, d aia2 5ip, d ai —a n Sf. After insertion of Eqs. (B9) and (BIO) into the Eq. (B8) and after their rearrangement we 
obtain as an intermediate result 



S£=d„ 



dC 



d(d ll ip r ) 



- d ai 



dC 



d(dp ai <p r ) 



d 



dC 



d(dfj, aia2 fr) 



dC 



d{d^ ai fr) 



d ai 



dC 



+ H n d ai ... an 



(-) n d ai - 



dC 



d(dfj, ai ... an fr) 



dC 



d^d^a^fr) 1 " d(d liaiai ... OCn fr) 

®£ _|_ _|_ / \ng ®£ 

d{d[ 1( j l (j 2 ipr) d{dfx ct 1 a 2 a± ■ ■ -a n fr ) 



Sifr (Bll) 
d ai Sf r 
d ai a 2 Sfr 



+ 



dC 



d ai - 



dC 



d{d^ ai ---a n fr) 1 d{d llai -a n a 1 fr) 



+ ■ ■ ■ + (-) n d ai . 



dC 



(T\ ■ --On OL\ ■ ■■OL n fr) 



As next step we replace 5p r by the total variation, Sxf, and insert the resulting equation into the total variation for the Lagrangian 
in Eq. (B6). Furthermore, we use 



d a CSx a = d^ (g£CSx a ) 



(B12) 



where g^, is a metric tensor. Eq. (B12) obviously holds, when the infinitesimal transformation for the coordinates, Sx^, concerns 
a constant displacement of the 4-vector x M . In case of rotations, where 8x^ depends on the coordinate x^ its self, Eq. (B12) still 
applies due to the antisymmetry of the tensor . These steps give us the final and general expression for the Noether theorem 



with respect to variations of the different fields and their coordinates 



5C=d u 



dC „ dC 



did^r) 1 d{d lloll ip r ) 
x (Sr'fir — d a ip r Sx a ) 

dC 

d{d l _ llJ1 Lp r ) 

x d ai (5 T ip r - d a ip r 8x a ) 



daia.2 



dC 



dC 



diduviaiVr) 



dC 



<9(<9 pCTlCT2 </5 r ) 



x d aitT2 (S T ip r - d a ip r Sx a ) 



+ (-) n d ai ... an 



dC 



d(d liai ... an tp r ) 



+ ■ ■ ■ + (~) n d ai . 



l ai ■■■ct-n 



+ ■ ■ ■ + (-) n d ai ... an 



DC 

d{^d^i G\(7<2 Q£l •••Qf n ^Pr ) _ 
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dC 



dC 



d(di ia - 1 ... a - n ip r ) 1 d{d IJirTl ... !Tn ai (f r ) 
x d ai ... an (S T tp r - d a ip r Sx a ) 

- g» a £Sx a \ . 



+ ■ ■ ■ + (-) n d ai ... am 



dC 



(7\---a n Qi ■■■a. n (fir) 



(B13) 



We consider now global phase transformations (e <C 1) 



Sx^ = , <p' r (x') = e le ip r (x) =>■ SlfirT = S(fi r — —ieifir 



(B14) 



and obtain the following relations for global phase transformations 



Sip r = — ieip r , 
d ai 5(p r = - ied ai ip r , 



Jipr = - ied aia2 ... an <p r . 



(B15) 



Therefore, the invariance of the Lagrangian under global phase transformations results to the continuity equation d^J^ = with 
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the current given by 
8L 



J" 



dC 



- Q a |- Q 

d(d^ r ) ai d{d^ ai (f r ) 
8C 



8L 

12 d(d^ aia2 ip r ) 



8 a 



8C 



+ 



d(dp ai <p r ) 8(8^ ip r ) 

8C 



(—) n 8 



dC 



(B16) 



+ --- + (-) n d ai ... ar , 



+ --- + (-) n d ai ... ar , 



8C 



<9£ 



8 ai <r 2 <Pr 



9£ 



8C 



d{d IJ , tT1 ...< Jn ipr) 1 9(3, 



CTi •••(J n a l 



tp r ) 



+ --- + {-) n d ai ... an 



dC 



8(8, 



p, <7\ ■■■(T n ai ■ ■ -a 



9 ai -a n 



Using Eqs. (10) for the tensors K,^ ai<J2 '" one arrives to the Expression in Eq. (9). 

The energy-momentum tensor is derived again with the help of Eq. (B13) for the case of constant 4-translations x' ^ 
x^ + 8x^. This means Sr'firix) = and following expression is obtained 



8C 



9 a 



8C 



8 a 



8C 



+ 



8C 



8 a - 

8C 



-■■■ + (-) n d al ... an 

+ --- + (-) n d ai ... an 
+ ■ ■ ■ + (-) n d ai ... an 



8C 



<Pr) 

8C 



fl <Tl(T2 ai 



d v ip r 



+ 



8C 



8C 



8C 



d(d llax ... an ip r ) 
-g^C^Sxv = 0. 
This leads to the continuity equation 



Fill \ "t" ' ' ' ( ) ^ ai ' an f)(fi \ 



with the energy-momentum tensor T^ v given by 



0„, — r + O n 



<9 M = 



8C 



d(d^ r ) " ai d(d^ ai ip r ) 01102 8(8^ <Pr) 



8C 



8C 



8(dfj.a iai ^r) 



8C 



(B17) 



(B18) 



+ {-) n 8 ai ... a „ 



8C 



+ ■■■ + (-) n 8 



8C 



+ --- + n n 8 ai ... a 



" 8(8^ ai <p r ) 
8C 



<Pr) 



8"<p r 



(B19) 



8C 



8jC 



8{d liai ...< Jn (p r ) Q1 9(9 M(71 ... CTn 



(-) n d a 



8C 



8(8, 



-gTC. 
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Above expression for the energy-momentum tensor can be written also in a more compact form resulting to Eq. (12). 



Appendix C: Preliminaries for the NLD formalism 



For the derivation of the Dirac equation, the current and the energy-momentum tensor in the NLD model we need the deriva- 
tives of the NLD Lagrangian with respect to the spinor fields and their higher-order derivatives. We evaluate them here up to 
second order and provide all terms of infinite series in the final expressions. Since we are interesting on the derivatives with 
respect to the spinor fields only, we start with the NLD Lagrangian without the meson-field contributions. That is 



9p 
2 



95_ 

2 



(CI) 



The application of the various higher-order partial derivatives with respect to the spinor fields and ^ to the Lagrangian density 
in Eq. (CI) proceeds with the help of the multiple Taylor expansions, see Eqs. (26). It is convenient to rearrange these series 
in ascending order with respect to the partial derivatives. With £ j = — Cf a and £ j = id aCf ti = 1, 2, 3, 4) where 
Cf = v> i I ^ one obtains for the expansion up to order n 



with the condition m + ■ ■ • + n 4 = n and 



j(n) 



: 



7(1) 



5 



-X>| 



gn 



(C2) 
(C3) 

(C4) 
(C5) 

(C6) 



The pairs between Latin and between Greek indices in above equations denote the summation over the multiple variables £j (i = 
1, 2, 3, 4) and over the 4-coordinates, respectively. In order to simplify the derivations in the following appendices, we skip the 
summation over the multiple variables. 

For the partial derivative of the Lagrangian density with respect to the spinor field $ only the zero-order terms in Eqs. (C2,C3) 
contribute. Therefore we obtain in detail 



9£ 1 1 

1 ? 
+ 2 9s6 



i 



1 



;9 P P f 



7m' 



7m' 



(C7) 



and similar for the first-order derivative with respect to the spinor field 



dC 1 -, U T t 1 



*^7 p f + *7 M T(i (0) 



(C8) 



Concerning the partial derivatives with respect to the first-order spinor fields d a ^> and d a ^ only the first-order derivative terms 
in Eqs. (C2,C3) contribute and we get 



dC 



=r = - \l ai iip + lg^d (1) i( ai - \gud {1) < Ql 7^a;„ - i.9 P d (1) < Ql 7 M f #p „ + < Ql r f<5 , (C9) 



2 
1 .- 



5 P c 

i 9(T a*zC Ql + \ 9 uio^r iC 1 d {1) + \ 9p p^t^ iC, ai - ^gsSr^iC 1 . (CIO) 
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In a similar way as above one evaluates the derivatives of the Lagrangian density with respect to the second-order partial 
derivatives of the Dirac spinors. In this case obviously the second-order derivative terms in Eqs. (C2,C3) are of relevance, and 
the result reads 



dC 



+^ (2) ^ri(" 2 ?«, (Cll) 

4 5ct a * 2T d(2) " 2^ < Q1 <" 2 d(2) ~ y | < Q1 <" 2 rf (2) 

+ ^.9^>^< ai < a2 rf (2) - (C12) 



With the intermediate results of this appendix we can now derive the relevant equations of the NLD model, i.e., the Dirac- 
equation for the spinor field "J in Appendix D as well as the conserved Noether current and the energy-momentum tensor in 
Appendix E. Furthermore, we will perform these derivations up to second order in the higher-order fields and for the isoscalar 
meson-nucleon interaction Lagrangians only, since the evaluation of higher-order terms and of the other meson-nucleon vertices 
proceeds in a similar way. We will insert then the remaining terms, i.e., the higher-order derivatives as well as all meson-nucleon 
contributions in the final expressions. 

In the following the terms containing the derivatives of the meson fields are not shown, since, they do not contribute on the 
mean-field level. 



Appendix D: Derivation of the Dirac-equation in the NLD model 



For the derivation of the Dirac-equation we start with the Euler-Lagrange equations of motion, Eq. (7), which read as 



dC 

dip r 



dC 



i=l 



d(d ai ... ai <p r 



(Dl) 



Up to second order in the partial derivatives of the spinor field they reduce to 

dC 



<9*~ ai d(d ai V) 



d(d ai d a . 2 V) 



(D2) 



We insert the various partial field derivatives, Eqs. (C7), (C9) and (Cll), into the second order Euler-Lagrange equations, 
Eq. (D2), and obtain 



+\g* o i C Q1 C 2 it ai it aa * - hg u d (2) ^ - C 1 C 2 i^ ai it aa * = o . 

We rewrite Eq. (D3) such to separate the series contributions from the standard terms and obtain following expression 



(D3) 



-m+ -g a aT$-- 



7^ 



73 



(D4) 



g<j o 
gujUf, 



d co) _ d d) c*i$ ai + 1 d& cn1) ai C 2 i$ a 



2! 



d (o) _ d (D £* H t ai + i d& C 2 i1> a2 7^} * = 
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In fact, if one would perform above procedure for all higher- order terms, one would obtain 

2* 



7^9'' -m+ -g a aV - -g u 7^ V 



(D5) 



1 



d (o) _ d (i) c a n1) ai + 1 ^ c Q1 ^ Q1 C 2 i^ a2 + ■■■ + (-)" 3 d (n) C Q1 ^ ai ■ • • C Q "^ 



2-9" w m 



2! n! 



7^ ^ * = . 



The infinite series inside the brackets in Eq. (D5) add together with the non-linear terms in the first line of Eq. (D5) for each 
meson-nucleon vertex. One arrives to the following Dirac equation for the spinor field * in the NLD model 



7m id" ~ g u 7^m^ - 



— m + q a a ~i> + qs tS T$ 



* = 0. 



(D6) 



This is the desired result, Eq. (28). Again the terms containing the derivatives of meson fields are not show in the above equation, 
since, they will not contribute in the final RMF expressions. 



Appendix E: Derivation of the Noether current in the NLD model 



As in the derivation of the Dirac equation, we consider in the following only terms up to second order in the field derivatives, 
and use Eqs. (C9), (CIO), (CI 1) and (C12). We start from the general expression in Eq. (9) for the nucleonic degrees of freedom 
which up to second order takes the following form 



dC n dC 
-dfi- 



L9(9 M *) p d(d^) 



* - * 



dC _ g dC 
9(9/,*) _ d{d^) 



dC 



9(9^*) 



9/3* - 9/3* 



dC 



ld(d^) 



(El) 



We rewrite now Eq. (El) by separating the terms between the different orders in the partial derivatives (the order is labeled by a 
subscript) 



j M = of + of 



The first-order contribution to Eq. (E2) reads as 

= -i 



dC T - dC 



9(9^*) 9(9^*)/ ' 
while the second-order contribution to Eq. (E2) takes the following form 

dC 



0(2) 



„ dC - dC 



d(dfdp*) 



d{d»d^)J \d{d^d^) 



9/3* - 9,3* 



dC 



9(9^9/3*; 



(E2) 



(E3) 



(E4) 



For the evaluation of the first-order contribution of, Eq. (E3), we insert Eq. (CIO) and its adjoint form, Eq. (C9), into Eq. (E3), 
and obtain 



Of = * 7/1 * - i g„ [a* C M * + * Cm * o 



(2) 

The derivation of the second-order contribution , Eq. (E4), proceeds in a similar way. We get 

Of =-\g.o^i% l? \c^n^^^\g^d^\c^id fi *a 



1 

2 9u 



^*7 q Cm d {1) * + * ^ 7 q *^c 



(E5) 



2! 



2! 



\ gu us * i dp^ Cm C P I s d {2) * - \ ft, * d (2) ^ Cm C^ 7 S i dp * 
i 9u] u s * 7 5 I C, Q P d (2) i dp * + \ g u * rf (2) i 9/3 ^ Cm C P 1* * 



W<5 

LOS 



(E6) 
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For each isoscalar meson-nucleon interaction we obtain now 4 terms, which differ between each other in the direction where the 
partial derivative operators act. We arrive to the following expression 



- i 9ui us W § Cm C 3 d (2) < 0^ * + \ 9u * rf (2) i 3/j | Cm C I s * w« 



(E7) 



The procedure is similar for the remaining higher-order derivative contributions. The evaluation procedure according Eqs. (E6) 
and (E7) for the third-order contribution, , would result in three terms for each vertex, for the fourth- order term, of, in 
four terms for each vertex, and so forth. 

Therefore, the resummation of all higher-order terms according 



J M = of + of + of + --- + o^ 



(E8) 



leads to infinite series for each meson-nucleon interaction. For instance, for the scalar-isoscalar meson-nucleon interaction we 
get 



i_ d (D ^ + ^ d (2) ^ (ai ^ + . . . + n d (n) fa CH <g c 2 --- iK n _ t c- 1 c 
1! 2! n\ 

~7T C M + I C M d {2) C 1 i$ ai +■■■ + (-)" -. e d (n) C 1 i$ ai C 2 i^« a ■ ■ ■ C- 1 i^a^ 
1! 2! n! 



* . (E9) 



Note that the n-th term, Oj 1 , contains (n — 1) partial derivatives and that it appears n-times. These terms can be also resummed 
to infinite series. Indeed, by considering, for instance, the Taylor-expansion of the operator V 



V=l- d« (a it" + ^ rf (2) Ca it" C/3 iti 

we obtain for the derivative of 2^ with respect to the operator- like argument 



. 



(E10) 



(Ell) 



In a similar way we obtain 7> = (7?/*)t. As in the case of the scalar operators T$ and i D, the new vector operators 7^ 
and act from the right side to the spinor field ^ and from the left side to , 3>, respectively. Below we will see that in the 
RMF approximation to nuclear matter both vector operators will give the derivative of the scalar operator V with respect to the 
single-particle 4-momentum p M . 

Collecting now all the contributions from Eqs. (E5) and (E7) under consideration of Eq. (El 1) we arrive to compact forms, 
e.g., the scalar-isoscalar meson-nucleon vertex from Eq. (E9) can be resummed as follows 



1! 



— ^ d (l) + ^ ^ d (2) ri + . . . + » ^ d («) ^ ."^ C a 2 . . . C a B _ 1 ^ 

1! 2! n! 



(E12) 



The evaluation method for the isovector channels of the NLD interaction Lagrangian proceeds in the same way. In total, including 
all degrees of freedom we obtain following compact expression for the conserved baryon current within the NLD formalism 



i 

+- 2 9 P 



" la - cr* 7^*1 + - q u , h t?^ 7 Q *a; Q - w a *7 a 7^* 



2 9u 



2 55 



(E13) 



which obeys the continuity equation 9„ J p = 0. 

Now we apply the RMF approximation to the general expression for the Noether current, Eq. (El 3). All bilinear operator-like 
terms are replaced by their expectation values relative to the nuclear matter ground state. Furthermore, we use for the spinor field 
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the plane-wave ansatz, Eq. (38), in order to evaluate the operators 2>, V, f> and J? M . Taking also into account the equations 
i~$t l \fj = p^xp and "J/i d M = — "fj/ 1 the current J M , Eq. (E13), takes following form in the RMF approximation 

= + g a (d»V)V)a - g u (0£Z>) 7 Q tf )w Q - <? p (¥ (c£2?) 7 Q f *)p Q + g s (* {d^V)f^)5 , (E14) 

with 9^ = g^-. This is the desired result, see Eq. (47). 

The derivation of the energy-momentum tensor proceeds in the same way as for the current, therefore, we skip further deriva- 
tions. 
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